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ABSTRACT 



' Aims. We investigate the cosmological growth of dark halos and follow the consequences of coeval growth for the accretion history of 
associated supermassive black holes. 

^SJ ■ Methods. The Press-Schechter approximation is used to obtain an analytic expression for the mean rate of growth of dark matter halos. Dark 

^ ' halo accretion rates are compared with numerical work and the consequences for understanding AGN evolution are described. 
C^l Results. The mean accretion rate onto dark matter halos is shown to have a simple analytic form that agrees with previous numerical work and 
■ that may easily be calculated for a wide range of halo mass, redshift and cosmological parameters. The result offers a significant improvement 
over published fitting formulae deduced from merger trees. We then consider the growth of associated supermassive black holes, and make a 
basic test of the simple hypothesis of "Pure Coeval Evolution" (PCE) in which, on average, black hole growth tracks dark halo growth. We 
demonstrate that both the absolute value of the integrated AGN bolometric luminosity density and its cosmological evolution derived from 
| hard X-ray surveys are well-reproduced by PCE. Excellent agreement is found at z i 0.5, although the observed luminosity density drops by a 

■ factor 2 compared with PCE by z = 0: black hole growth appears to decouple from halo growth at low redshifts, and this may be related to the 

' phenomenon of "cosmic downsizing". Overall, AGN evolution appears either to be caused by or to be closely linked to the slow-down in the 
growth of cosmic structure. We also discuss the mean Eddington ratio averaged over all galaxies, which is predicted to show strong evolution 
Q ■ to higher values with redshift. 
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1. Introduction Yet this picture is not easily in accord with modern ideas 

of galaxy formation, since if the matter content of the universe 
is dominated by cold dark matter we expect galaxies to grow 
One of the great mysteries in the study of active galaxies and hierarchically, and we expect the mass function of black holes 
QSOs is the physical origin of the strong cosmological evo- at the centres of galaxies to increase with time - what then is 
lution in their space density at a given luminosity. For many the mechanism that allows the mass function to increase but 
years it has appeared that, at least in broad terms, the evo- the luminosity to decrease with cosmic epoch? One obvious 
lution was best described as "pure luminosity evolution" in possibility is that black holes may increase in mass but with 
which QSOs appeared on average to h ave dimmed with time a decreasing mean accretion rate. This picture has been incor- 
rather than chang i ng in space density llMarshall et all fl983t porated into models such as those of |Haehnelt & Rees| \l99j) 
iBovle et all \l9M ICroom etall Eooi . This is in accord with and|Kauffmann & Haehnelt| (|2000|) in which simple parameter- 
recent evidence that almost all galaxies with a massive spheroid isations of such cosmological evolution in accretion rate were 
compon ent contain a supermass i ve black hole at the present allowed, and it was shown that with a suitable choice of pa- 
epoch dMagorrian et all fl998t iFerrarese & MerritA EoOp] rameters the evolution in the QSO optical luminosity function 
Gebhardtetal., 2000; Tremaine et al., 2002; Onkenetal., could be matched. The physical origin of this variation in ac- 



2004), the great majority of which must be largely inactive to- cretion rate has not y et been determined, however, and expla- 

day. If these black holes were luminous in a phase of accretion nations ran 8 e from a systematic depletion in available gaseous 

and growth at higher redshift then their luminosities must in- material with cosmic epoch to long-term variations in the ac- 

deed have declined with time cretion process itself: these latter explanations would of course 
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only work if black holes were formed at high redshift and the 
luminosity evolution that has been observed were a reflection 
of the variation in accretion history of individual black holes. 

More recently however it has become clear that optically 
luminous QSOs are a short-lived phenomenon, at least com- 
pared with the Hubble time. Measurement of t he clustering 
of QS Os in the 2dF QSO Redshift Survey (2QZ: ICroom et all 
2005) has shown that their clustering amplitude does not in- 
crease with time as would be expected if QSOs were long-lived, 
given that cosmic structure is expected still to be growing at 
z ~ 2 in a Universe with a cosmological matter density param- 
eter Q.m ~ 0.3. This implies that we must look for a universal 
variation in accretion rate rather than anything intrinsic to an 
individual accretion "event". 

The simple picture of QSO luminosity evolution itself has 
now been shown to be more complicated when active galaxies 
are selected at X-ray wavelengths. First, it now seems likely 
that the hard X-ray background is produced by X-rays emitted 
from active galaxies, but in that case a large fraction of those 
active galaxies must be highly absorbed, with equivalent X-ray 
absorp tion co l umn densities in neutral hydrogen > 10 23 cm' 2 
dRisaliti et all Il999t IComastri et all EoOlt llJeda et all 120031) . 
These active galaxies would be likely to be optically obscured 
also. A second significant discovery from X-ray surveys how- 
ever has been that low luminosity X-ray selected active galax- 
ies, that would be classified as Seyfert galaxies were they op- 
tically selected, display cosmological evolution in their space 
density that appears to have a maximum at l ower red shift 
than their high luminosity QSO counterparts ( Steffen et al. , 
2003HCowie et alll2003UlJeda et allfcoosllzheng et alll2004 : 
Barg eret all 12005). 't may be to some extent that the identi- 
fication of these weak active galaxies may be less complete at 
higher redshifts, but early indications are that there is a signifi- 
cant shift of the space density maximum to lower redshifts with 
decreasing luminosity. 

Hence it appears that there is much still to understand both 
observationally and theoretically about the evolution of active 
galaxies. Nonetheless, the basic inference that there must have 
been evolution in accretion rate seems inescapable. In this pa- 
per we test the extent to which the growth of black holes may 
be linked to the growth of galaxies, and argue that cosmolog- 
ical evolution in accretion rate can be understood primarily as 
arising from the cosmological evolution in the rate of accretion 
of matter onto galaxies. In section 3 we present a new calcu- 
lation of the accretion rate onto dark halos using the extended 
Press-Schechter approximation. We find that the analytic ex- 
pression agrees well with both our own and previous numerical 
estimation of the accretion rate calculated from merger trees. 
The analytic approach here circumvents the need for such nu- 
merical estimation, and we also find that it provides mathemat- 
ically better behaviour than fitting functions to the numerical 
results that have previously been proposed in the literature. 

In section 4 we argue that coeval growth of black holes with 
their associated galaxies and dark halos leads to a simple ex- 
pression for the total mass accretion rate onto all black holes. 
We calculate the expected integrated luminosity density arising 
from accretion onto black holes and show that it agrees remark- 
ably well with observation. When expressed as the Eddington 



ratio (the ratio of the actual mass accretion rate to the rate re- 
quired to attain the Eddington limiting luminosity) it shows 
significant cosmological evolution almost independent of halo 
mass and without depending strongly on the choice of cosmo- 
logical parameters. 

2. Coeval growth of supermassive black holes and 
their host galaxies and halos 

One hypothesis, whose consequences will be followed in this 
paper, is that the population of supermassive black holes, that 
in the nearby universe inhabit the nuclei of massive galaxies, 
grew coevally with their host galaxies and associated dark ha- 
los. Before launching into detailed calculation, in this section 
we first consider the rationale for testing this hypothesis. 

First, given present data on the black hole/bulge relation- 
ship, some degree of coeval growth seems hard to avoid, at least 
for galaxies with measurable bulge components. In all massive 
galaxies where both a bulge mass and a black hole mass have 
been measured dynamically, ther e is an extremely tight corre- 
lation between the two jpebhardt et all l200ol iTremaine et all 
l2002HOnkenetall f2004). At lower masses, with velocity dis- 
persion cr < 100 kms -1 , the correlation is less well-established, 
but does appear to extend at least down to cr ~ 30kms _1 
(Bar th et all 12005 s ) . When the bulge mass is determined from 
velocity dispersion measurements, the cosmic scatter in the 
relation appears smaller than the me asurement erro r s, and 
probably is less t han a factor three jTremaine et all |2002; 
IOnkenetall E004). But it is thought that galaxies have grown 
hierarchically over a wide range of cosmic epoch, and with 
galaxies at any observed redshift having accumulated their 
masses at many different times. In order for the tight present- 
day black-hole/bulge relation to exist, the process of formation 
of black hole and galaxy must have been interrelated, the sim- 
plest explanation being that they grew in mass at the same time. 

Second, the luminosity emitted by AGN is thought to 
arise from accretion, and therefore is a signatur e of o ngo- 
ing black hole growth at low redshifts. lUeda et all 12003) and 
Marconi et al. (20Q4J) have shown that both the X-ray back- 
ground and the local mass density in black holes are consis- 
tent with being created by the luminous phase of accretion 
visible in hard X-ray AGN surveys, and that the bulk of that 
black hole growth occurs at low redshifts (z < 3). The pic- 
ture of supermassive black holes forming at high redshifts and 
then remaining largely unchanged since then is not consistent 
with the observed AGN luminosity density and the inferred lo- 
cal black hole mass density, unless the radiative efficiency of 
the AGN luminosity we can see is unfeasibly high (e >> 0.1, 
iMarconi et alJl2004l) . So it appears that black holes have been 
continuing to form during the cosmic epochs in which dark ha- 
los and their associated galaxies have also continued to grow. 

Previous AGN models (e.g. lKauffmann & Haehnellll200oT) 
have assumed a relationship between dark halo mass and black 
hole mass. It is hard to see how the black hole could know what 
mass of halo it is in unless there has been some causal link 
between halo growth and black hole growth. The link does not 
need to be direct, and it may be that feedback processes have 
an important role in regulating the black hole/halo relationship, 
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but in effect, these models have implicitly assumed that coeval 
growth has occurred. Whether coeval growth is still occurring 
in the present-day universe is another question, and one that we 
shall attempt to investigate in this paper. 

The growth of galaxies by hierarchical mergers implies that 
there is a contribution of mergers to black hole growth, too. 
Thus we expect the mass function of black holes to be de- 
termined by hierarchical merging, as is the mass function of 
galaxies. Mergers of black holes cannot change the integrated 
irreducible mass in black holes, however, so the integrated 
mass in black holes should depend on the total amount of mat- 
ter accreted and not on the merger history. The importan ce of 
mergers does not violate the argument originally due to Soltan 
( 1982) that black holes in the local universe are produced by 
accretion at earlier cosmic epochs. The growth of supermas- 
sive black holes through mergers does in principle allow high 
black hole masses to be attained at rates faster than the Salpeter 
rate without violating the Eddington limit, which may help us 
to understand how supermassive b lack holes can exist at high 
redshift, z ~ 6 (e.g. lWillott et all2003l) . 

In practice we might expect the phases of significant mass 
accretion in the life of a black hole to be correlated with periods 
of galaxy merging: it has long been suggested that mergers may 
drive matter in to galaxy centres and trigge r phases of black hole 
accretion (e.g. Barnes & Hernquist 1992). This picture creates 
an attractive link between galaxy mergers and prompt AGN 
activation but is not required in what follows. 



3. The mean accretion rate onto dark matter halos 

3.1. Structure growth from the Press-Schechter 
approximation 

If the matter content of the universe is dominated by cold dark 
matter, then dark matter h alos a ssociated with galaxies form 
hierarchically. iBond et al] dl99 l|) have shown tha t a mo re rig- 
orous treatment of the work of lPress & Schechterl Jl974h ("ex- 
tended Press-Schechter", hereafter EPS) can be used to ob- 
tain information about the build-up of structure, and the result- 
ing evolving mass functions agree with the results of N-body 
simulations, especi ally if some additional modification of EPS 
theory is allowed dSheth & Tormeii 1 19991 ISheth et all 1200 it 
Sheth & Torrrienl l20Q2h . 

In this picture of hierarchical galaxy formation, at any 
moment in time, any given overdensity is increasing in mass 
through the process of accretion of matter. There has been 
much discussion in the literature about whether one can use 
EPS theory to further analyse thi s accr e tion in terms of 
merger events (Lacev & ColeL Il993l 1 19941 ICole et all [2000; 
iBenson et all 120051) but we are not interested here in attempt- 
ing suc h a detailed view of the build-up of galaxies. In previous 
papers JPercival & Millerl 1 19991 IPercival et all l2000l) the for- 
mation rate of halos within EPS theory has been calculated, and 
hence if the space density of halos is also known this may be 
converted into a mean accretion rate onto halos as a function 
of their mass and redshift, measuring the overall accretion of 
matter, not broken down into individual events. 



In those papers it was shown that a joint distribution func- 
tion in mass M and cosmic time t of new overdensities, F(M, t), 
could be defined that is related to the mass function, F{M\t), by 
the relation 

dd c 



F(M,t)dMdt = F(M\t) 



dt 



dMdt, 



(1) 



where 6 c (z) is the linear- theory overdensity at redshift z — 
that would be required for an overdensity to have collapsed at 
some redshift z > 0. For an Einstein-deSitter universe 



dt 



1.686// (1 +z) 



2.5 



(2) 



In Appendix A we discuss alternative methods of calculating 
\d5 c /dt\ for other choices of cosmological parameters, and de- 
scribe an analytic approximation based on the linear growth 
factor of cosmological perturbations. 

However, the normalisation of the joint distribution F(M, t) 
is problematic, as discussed by IPercival et all (120001). In 
Section l3~2l we follow a different approach and calculate the 
mean accretion rate analytically from the EPS conditional mass 
function. This is complemented in Section l3~3l bv a numerical 
derivation of the rate from Monte-Carlo realisations of the EPS 
process. 

3.2. Analytic derivation of the mean halo accretion rate 

Lacev &Colel il 19931 1 1994) give the conditional probability 
that a halo of mass M at cosmic time t was previously a halo of 
lower mass M' at earlier time t': 

— (M',t'\M,t) = 
dM' K ' ; 



A6 C 



VMAo- 2 ) 



.2^/2 



do- 2 



dM' 



exp 



(A<5 e ) 2 



2Acr 2 



(3) 



where A6 C = 5 c (f) - 6 c (t), 5 c (t') > 5 c (t), and where A<x 2 
it'' - cr 2 with cr' 2 being the variance on the mass scale M', cr 2 
the variance on the mass scale M and cr' 2 > cr 2 for M' < M. 

Hence the expectation value for the halo's increase in mass 
AM = M - M' is 



{AM) 



(M - M')A6 C 



do 2 



dM' 



f 

Jo VMAo- 2 ) 3/2 

Making the substitution y 2 
give 

{AM) = 



exp 



dM'. 



(4) 



2(Acr 2 ) 

1/Acr 2 this can be rearranged to 




MA5 C \ l — 



M'(y) 
M 



expf-y 2 (A<5 C ) 2 I2]dy. 



(5) 



As y -> oo, M'(y) -> M and (l - ^jp-) — > 0. Hence we can 
make A5 C sufficiently small that the exponential term can be 
ignored in this integral, and in the limit At = (t - t') — > 0, 
A6 C — > and 



(AM) MAS C \j — 




M'(y) 



dy 



M- 



V2^ 



M 

M'\ dAo- 2 
M I (Act 2 ) 



2YV2- 



(6) 
(7) 
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where Aa 2 — > oo as M' — > for CDM power spectra and 
Act 2 — > as M' — > M. The kernel of this integral was de- 
rived bv lCole et alJ J2000h but those authors were interested in 
considering the distribution of merger events rather than inte- 
grating to find the overall mean accretion rate, which forms the 
basis of this paper. 

Hence we can write the mean rate of mass accretion as 



(«) = 

where 

/ 



,. (AM) 
hm 

A5 c ->0 A5 C 



d5 c 
dt 



M 



dt 



f(M), 



1 f 00 / M' 



dAo- 2 



(Acr 2 ) 3/2 ' 
Integrating by parts we find 

M'=M 



(8) 



(9) 



f(M) 



KM 



1/2 



M'=0 




m- 1/2 — . 

V / M 



(10) 



The first term is zero provided 



/(M)= vf^r^ i/2 ' M '- 



> 0, leaving 



It should be noted that this particular definition of mass ac- 
cretion rate is not the only one possible. One could also cal- 
culate the accretion rate by considering the future accretion of 
halos and taking the limit of small time intervals: this produces 
a different analytic form, involving an integration over masses 
larger then the mass being considered, and highlights that the 
mass accretion process for a mass M at time t is discontinuous 
within EPS theory. In section l3~3l we also consider the relation- 
ship of this "instantaneous" mass accretion rate to the process 
of the build-up of halos on longer timescales: the mean accre- 
tion rate is a function of the time interval being considered, as 
may be seen from inspection of Eq.[5] 

The function f{M) may be evaluated numerically for an 
assumed matter power spectrum. For cold dark matter power 
spectra it depends only weakly on mass and varies little with 
shape of the power spectrum. It varies inversely with the nor- 
malisation parameter erg. Fig. [2 shows the variation of f{M) 
w ith halo mass for a AC DM power spectrum approximated as 
in Efstath iou et al.l i 19921) for a range of values of the shape pa- 
rameter F = Om/i and with fixed cr 8 = 0.74, h = 0.7. The rela- 
tionship between variance <x 2 and mass was calculated assum- 
ing a spherical top-hat smoothing function, although we should 
note that strictly Eq. [3] is derived assuming smoothing with a 
top-hat function in fc-space. At a halo mass of 10 l 2 5 h~ l M n , ap- 
propr iate for the halo masses of luminous QSOs JCroom et all 
2005]), f(M) varies only by a factor 2 over the range 0.1 < 
Qm < 0.4 for h = 0.7. The variation with mass is also weak 
for masses in the range appropriate for massive galaxies: at 
Qm = 0.3, h = 0.7, f{M) varies by a factor 2 over the range 
10<log 10 (M/M G )< 13. 
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Fig. 1. The variation of the function f(M) with logio(halo 
mass) as a function of the power-spectrum shape parameter T, 
with h = 0.7 and cr 8 = 0.74. Curves from top to bottom have 
T = 0.07,0.14,0.21,0.28,0.35,0.7; solid curves show the cal- 
culated function, dashed curves show the fitting function, equa- 
tjorJPH 

For convenience we provide a fitting formula for f(M). 
Over the range of F and mass shown in Fig.^ we adopt the 
relation 



(11) log 10 (/[M])=A(r) + (fi-A(F)) 



lo gio 



(/i/0.7) 2 M^ 
10 17 M o ; 



(12) 



with B = 1.117, B = 4.82 and A (T) = -0.7349 - 
0.9808(log I0 F)-0.2055(log 10 T) 2 . This function fits to an r.m.s. 
accuracy of 1.6 percent over the range 0.07 < F < 0.7, or 
0.1 < Oo < 1 for h = 0.7, and is shown by dashed lines in 
Fig-E For cr 8 = 0.84 the best fitting function parameters are 
B = 1.072, B = 4.57 and A (F) = -0.8665 - 1.002(log 10 F) - 
0.2093(log 1() F) 2 . This formula may be used together with the 
analytic approximation to \dS c /dt\ described in Appendix A to 
calculate the mean cosmological accretion rate for a wide range 
of cosmological parameters. It is worth noting that the accre- 
tion rate has a simple dependence on mass and redshift, and 
either weak or at least linear dependence on cosmological pa- 
rameters. 

3.3. Numerical calculation of the mean halo accretion 
rate 



To validate the derivation presented in Section 13721 we now 
compare the analytic instantaneous rate with the rate calculated 
from Monte-Carlo realisations of the EPS process. In EPS the- 
ory, the overdensity recovered if a realisation of the density 
field is smoothed on different scales around a fixed point can 
be matched to a one-dimensional Brownian random walk. The 
collapsed mass at a given time is given by the trajectory's first 
upcrossing of a barrier at the critical overdensity corresponding 
to the chosen time. Because we can only calculate finite repre- 
sentations of these trajectories, there are a number of numerical 
issues that complicate their use. In particular 

- the discrete nature of the numerical trajectories means that 
we can miss upcrossings of a barrier between steps in the 
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redshift 



Fig. 2. Comparison of the mean mass accretion rate calculated 
analytically (solid curves) and measured from Monte-Carlo re- 
alisations (solid symbols), shown as a function of redshift for a 
variety of masses. 



realisations. If the step size is chosen to be too large, then 
this effect becomes particularly apparent for the first step in 
each trajectory, where most upcrossings are expected. 
- the number of steps must be finite, so we will always miss 
upcrossings at high a 2 (small mass). 

By appropriate choice of the step size, we can minimize 
these problems. We have created a sample of 10 6 trajectories 
each consisting of 2 13 steps of varying size in cr 2 . A vary- 
ing step size was chosen to more closely match the expected 
upcrossing distribution and allow a larger range in cr 2 to be 
probed. The barrier height was set at 4 times the initial step 
size to minimise problems owing to missed upcrossing in the 
first few steps of the trajectories. Because the trajectories can 
be arbitrarily scaled in 6 2 /cr 2 , we can use this single set of tra- 
jectories (and associated upcrossings) to determine the mass 
accretion rate for any halo mass and epoch. For the masses and 
times chosen, we have determined that the problems owing to 
the finite extent of the trajectories are negligible. By comparing 
with the expected conditional mass function (Eq. |3, we have 
also optimised the step size to reduce problems owing to the 
discreteness of the realisations. 

Fig-El shows the result from the Monte-Carlo measure- 
ments, calculated from the average mass accreted within A6 C = 
0.032, corresponding to 0.0035 Gyr at z = 7 and 0.5 Gyr at 
z — 0. Decreasing A6 C further does not significantly change the 
recovered rate. The agreement with the analytic instantaneous 
rate is extremely close. Any differences may be revealed by 
comparing f(M) with the equivalent quantity from the Monte- 
Carlo realisations by factoring out the analytic linear depen- 
dence on mass and d6 c /dt (Fig[3J. The maximum difference 
between the analytic and numerical results is < 5% and is likely 
to be caused by remaining problems owing to the discrete na- 




10 12 14 

log 10 (mass / M ©) 

Fig. 3. Comparison of f(M) (solid curve) with that deduced 
from the Monte-Carlo realisations by factoring out the linear 
dependence on M and d6 c /dt. 

ture of the numerical trajectories. Note that the curves in Fig.|2] 
show the immediate mass accretion rate onto a halo which has 
a specified mass at a specified redshift, not how the mean ac- 
cretion rate varies with time for an individual halo. 

This comparison has shown that the results obtained from 
numerical realisations of the Press-Schechter process are in- 
deed consistent with the analytic results that we have obtained. 
This type of realisation forms the basis of the commonly-used 
"merger trees" and in the next section we compare our results 
with earlier numerical attempts to understand the growth of 
dark matter halos. 

3.4. Comparison with previous numerical calculations 

IVan den Boschl J2002I) has calculated the average mass accre- 
tion histories of dark matter haloes based on the extended 
Press-Schechter formalism and t he jV-branch me rger-tree con- 
struction algorithm o flSomerville & KolatdJ 19991) . He provides 
a fit to the numerically-obtained mass accretion histories by fit- 
ting a function of the form 



log 



/M(z)\ 


= -0.301 


log(l +z) 


\ Mo j 




log(l +Zf) 



(13) 



where M(z) is a halo's mass at redshift z, Mq is its present-day 
mass, and v and zj are non-independent parameters. The best- 
fit set of parameter values is v = 1.362 + 1.858 log(l + zj) — 
0.032 log [Mo/CIO'^-'Mq)]. z f is obtained from the root of 
equation: 



(14) 



S c (z f ) = SM + 0.477 ^2 [crHfiMo) - cr 2 (M )] 

where is set to a value f, - 0.254. l^n den Boschl J2002I) rec- 
ommends differentiating M(z) in order to find the mean mass 
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Fig. 4. Comparison of the mean accretion rate from this paper 
(dott ed lines) with the tw o alternative numerica l fits proposed 
bv IVan den Boschl d2002l) (solid lines) and bvh^chslgr^^lj 
(2002) with parameter values estimated by IVan den B osch 
(2002) (dot-dashed lines). See the text for more details of the 
calculation of these curves, which ar e shown for va l ues of 
log 10 (M(z)/M o ) = 8, 10, 12, 14. The IVan den Boschl ^002) 
curves becomes unphysical at low redshift and depart signif- 
icantly from the analytic results. Both numerical fits become 
unphysical at high redshift and are truncated outside the region 
of validity, as discussed in the text. 



accretion rate. As discussed by IVan den Boschl J2002), con- 
struction of such accretion histories numerically is not trivial, 
and depends on the definition of a halo's progenitor within a 
merger tree. However, it will be interesting to compare the re- 
sults of the analytic result presented here with these earlier nu- 
merical attempts; 

Wechsler et al 1 ll2002l) have also discussed the growth of ha- 
los in merger trees, and find a reasonable parameterisation of 
the accretion histories of halos to have the form: 



M(z) = M Q e 



(15) 



where a is a parameter. In fact, comparison with equa- 
tionOO suggests this is a good approximating function for 
I Wechsler et alJ J2002I) to have chosen. From equation|S]we find 



dz dz 



(16) 



In the case of an Einstein-de-Sitter universe dS c /dz = 1.686 
and a solution for M(z) of the form of equation[2] immedi- 
ately follows if we ignore the slight redshift dependence intro- 
duced by the f(M(z)) term. This solution will not be correct for 
other cos mologies, but should s till qualitatively be of the cor- 
rect form. lVan den Boschl (1200 2) has related the parameter a to 

the previous parameter zj by the relation a — \zf/ 1 .43 1 

In Fig|4]we compare our accretion rate with those obtained 
from these two alternative fits to the numerical mass accretion 
histories. The solid line shows the mean accretion rate obtained 



by differentiating eauationlT3*land the dot-dashed line from dif- 
ferentiating eauation[T5l bo th with the prescr iptions for the 
free parameters given by V an den Boschl 12002) that are quoted 
above. The dotted line shows the mass accretion rate from 
equation|8](note that, as in Fig.|2] each curve corresponds to a 
fixed value of mass at each redshift: i.e. these are curves of con- 
stant M(z) and not constant Mo). The comparison is made for 
several values of halo mass: log 10 (M(z)/ M ) = 8, 10, 12, 14. 
There is excellent agreement between the three curves over 
most of the range of redsh ifts. However, at the lowest redshifts 
the IVan den Boschl J2002I) accretion rate shows significant de- 
partures, and it is straightforward to see that this arises from 
the choice of fitting function in eauationlOl in the limit z — > 0, 
dM/dt — > if calculated from that equation, for parameter val- 
ues v > 1. Such behaviour is unphysical, it implies that all 
halos at the present day have stopped growing. This is purely 
an artefact of the choice of fitting function and leads to signif- 
icant departur es from the analytic ca lculation for z < 0.5. The 
result from the lWechsler et al.N2002h form consistently results 
in slightly higher values of accretion rate, but predicts the same 
behaviour as the analytic calculation at redshifts approaching 
zero. 

We also find a further limitation of the numerical pre- 
scriptions, which is that the calculated growth histories 
become unphysical at high masses and redshifts, with prob- 
lems arising when the present-day mass Mo > 1O 15 M . 
The problem is that the M(z) functions for differing masses 
cross at high redshift, such that a high mass halo might 
be predicted to have had a lower-mass progenitor than a 
lower mass halo. For example, the main progenitor mass 
of a M - 1O 16 M halo is predicted by eauationlP3l to 
be smaller than that of a M = 1O I4 M h alo at redshifts 
z_Z 2.5. We roughly estimate that the IVan den Boschl 
d2002l) fit can be used up to (log IO (Mo/ M ),z) = 
{(14.5, 1), (13.3, 2), (1 1.7, 4), (10.8, 6), (10.2, 8.)}, where 
the pairs of numbers denote the locus on the Mo,z plane 
beyond which the prescription becomes unreliable. The 
same problem is found when using the parameterisation of 
eauation[T5l with the range of valid redshift and present-day 
halo mass moved to lower values: (log IO (Mo/ M ),z) = 
{(13.3, 1), (1 1 .8, 2.), (9.8, 4), (8.3, 6), (7.7 , 7)). G iven the re- 
lationship between the IWechsler et all J2002) and analytic 
functions, it can be seen that the problem in this case arises 
from the prescription for the value of a rather than being an 
intrinsic problem with the choice of numerical function. 

Overall we conclude that the analytic and previous numer- 
ical results are in reasonable agreement in the regions of mass 
and redshift space where the numerical parameterisations are 
reliable, but that there exist significant regions of that space 
where the previous numerical fits become unusable. 

3.5. How good is extended Press-Schechter theory? 

So far, both the analytic calculation and the numerical reali- 
sations have relied on the EPS formalism. This formalism ef- 
fectively amounts to being a linearised approximation to the 
true process of non-linear growth of structure, and a priori we 
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should have no right to expect that this will be an accurate 
approximation. How can we be confident that the work pre- 
sented here has relevance to the real process of halo growth? 
There have been many comparisons of the calculation of the 
halo mass function predicted by EPS with that obtained from 
N-body dark matter s imulations an d ; with some caveats, good 
agreement is found. Lfenkins et alJ d200lh show the level of 
agreement when compared with Virgo consortium simulations. 
The mass functions agree to a factor 1 .6 over a wide range of 
halo mass. In this paper we investigate the rate of growth of 
dark halos: if the EPS growth rate were significantly wrong it 
would be hard to then reproduce a good estimate of the mass 
function, and it seems likely therefore that the calculations pre- 
sented here are correct when compared with dark-matter-only 
simulation s to this level of accuracy. 

ISheth & Tormenl 17999) have shown, however, that the EPS 
mass function may be modified to produce even closer agree- 
ment with the N-body s imulations. Th is work has been devel- 
oped further bv lsiieth et alJlEool andlSheth & Tormenl d2002) 
and confirmed in the Jenk ins et alJ |200 1 ) simulations. The 
principal modification is to reduce the critical density for col- 
lapse, S c , by a factor yfa where a is determined from compar- 
ison with the simulations to have a value a ^ 0.7. Insertion of 
this factor into equation|3]would simply imply that /(M) is also 
reduced by a factor -\fa ^ 0.84. This causes a slight shift in the 
curves plotted here but otherwise has no effect on the redshift 
or mass dependence. 

We note that the merger trees o flVan den Boschl J2002) have 
also been compared with N-body simulations, and reasonable 
agreement was found there too, albeit with a tendency for halos 
to form at higher redshifts which could also b e reproduced by 
introd uction of the factor ~\/a. In the case of the lVan den Boschl 
( 2002) simulations the best-fitting value for -\fa varied over the 
range 0.82 < a/a < 0.94, depending on halo mass. Hence 
it may be that there is some systematic and possibly mass- 
dependent departure of the EPS accretion rate from the "true" 
(dark-matter N-body) rate, but most likel y at a level < 20% 
in accretion rate. Finally, we note that Benson et al. ll2005l) 
have also argued that the distribution of mass discontinuities 
in EPS does not allow a self-consistent interpretation of these 
as merger events. This casts doubt on one of the fundamental 
assumptions of EPS-based merger trees, but should not affect 
the calculation of the mean accretion rate, which does not rely 
on that assumption. 



4. The accretion history of supermassive black 
holes 

4.1. Pure coeval evolution 

We now move on to consider what effect the evolution in mass 
accretion might have on the ultimate accretion onto a super- 
massive black hole lying at the heart of a massive galaxy. In 
principle, understanding not only the formation of galaxies 
within dark halos but also the accretion of baryonic matter onto 
a central black hole requires complex and poorly-constrained 
physics. However, we could make a very simple assumption 
for the net effect of these complex processes, and see how well 



the predictions compare with observation. The assumption we 
shall make is that, on average, the relative growth rates of cen- 
tral black holes track the growth rates of their associated dark 
halos (section|2li. That is, 



1 dM BH 
M B h dt 



1 dM H 
M H dt 



</(M H )> 



d& c 
dt 



(17) 



where Mh is a dark halo's mass, Mbh a black hole's mass, 
and where () denotes averaging over an ensemble of halos and 
black holes. Note also that we haven't needed to specify what 
fraction of mass ends up accreting onto a black hole - that cal- 
culation would indeed need some detailed physics - the hy- 
pothesis simply assumes that the ensemble-averaged relative 
growth rates are the same. We shall term this simple hypothe- 
sis "Pure Coeval Evolution" (PCE). 

A simple extension of this hypothesis would allow some 
non-linear evolution of the relation between black holes and 
dark halos. EQuationfTTl could be modified by a factor a such 
that 



1 c/Mbh 



Mi 



BH 



dt 



1 dM H 
M H dt 



a {f (Mm)) 



dS c 
dt 



(18) 



which w ould lead to a rel ation of the form Mbh k M^ as in- 
ferred bv lFerraresel 12002). We can see that the effect of such a 
non-linear term would be to modify the accretion rate calcula- 
tions only by a constant factor. 

4.2. The AGN luminosity density and its evolution 

In this section we calculate the observationally-determined 
AGN luminosity density and compare with the value predicted 
by PCE. We define the integrated luminosity density to be 
the luminosity emitted during the black hole accretion pro- 
cess per comoving cubic Mpc, summing over AGN of all lumi- 
nosities, as previou sly calculated by a number of authors (e.g. 
iMarconi etai1l2004 . The quantity predicted by the PCE hy- 
pothesis is the bolometric luminosity density, so we will need 
to apply bolometric corrections to observed quantities. 

4.2.1 . The bolometric luminosity density derived from 
X-ray surveys 

We shall first e s timate the bolometric luminosity density from 
the llJeda et all d2003h AGN luminosity function. It has re- 
cently been recognised that hard X-ray surveys currently pro- 
vide the most com plete way of selecting AGN of all types, and 
lUeda et al J ( 120 03) incorporate into their analysis correction for 
X-ray obscuration in the detected AGN (although if there is a 
population of extre mely absorbed, Com pton-thick, AGN these 
will be missing from Ueda et al. 2003's census - see below). 

Fig.|5] sho ws the bolometric luminosity density Pbol de- 
rived from the llJeda et all J2003h LDDE model at 2 - 1 ke V in- 
tegrated over the absorption-corrected luminosity range 10 40 < 
Lx < 10 48 erg s -1 . The X-ray luminosity density has been con- 
verted to a bolometric luminosity density and a small correc- 
tion for missing Compton-thick AGN has also been applied, 
as described below. The luminosity density is not sensitive to 
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the precise limits of integration: increasing the lower limit to 
Lx > 10 42 ergs~' decreases the z — luminosity density only 
by 0.08 dex without changing the maximum value (to 0.01 dex) 
or its redshift. Decreasing the upper limit to L x < 10 46 erg s _1 
has no effect to within 0.01 dex. Most of the bolometric lumi- 
nosity density is produced by AGN with Lx ~ 10 44 erg s" 1 . 

Also shown in Fig. [5] are estimates of the uncertainty in 
Pbol- To calculate the uncertainty , we fi rst refit the binned 
data shown in Fig. 1 1 o flUeda et al.ld2003l) . kindly provided by 
Y. Ue da, with the same six-parameter model (see lUeda et all 
120031 for details of the model), but with the normalisation pa- 
rameterised instead as the integrated luminosity density. We 
then estimate confidence intervals on the integrated luminos- 
ity density for each of the five redshift slices (0.015-0.2,0.2- 
0.4,0.4-0.8,0.8-1.6,1.6-3.0) in turn, marginalising over the re- 
maining five parameters of the fit. This procedure thus includes 
the full statistical uncertainty in the data and includes an ele- 
ment of uncertainty arising from the fact that to cover the full 
luminosity range requires some extrapolation of the model, al- 
though with the caveat that only the set of model functions 
that may be parameterised by the six-parameter function are 
allowed. The best-fit values for each redshift slice and the de- 
duced uncertainties are shown in Fig.|5] Some po ints deviate 
somew hat from the nominal curve deduced from lUeda et alJ 
(2003)'s best fit: this is likely because their best fit was obtained 
from a maximum likelihood fit to unbinned data whereas our 
points and errors were evaluated from fits to the binned data. 
We nonetheless expect that the size of the error bars should 
also be a good indication of the size of the 68 percent confi- 
dence region for the likelihood-determined function. 



4.2.2. X-ray absorption and the Compton-thick fraction 

It is recognised that 2 - lOkeV X-ray flux es must be cor- 
rected for absorption, and lUeda et all d2003h derive an esti- 
mate for the luminosity-dependent distribution of absorption 
column which has been incorporated into their analysis. In ad- 
dition, it is likely that some fraction of AGN are Compton- 
thick, that the Compton-thick fraction is a function of intrinsic 
luminosity and tha t even hard X-ray su rveys miss these objects 
JUeda et alll2003l) . IUeda et all J2003) assu me the (very uncer - 
tain) correction, based on the results of Risaliti et al I II 19991) . 
that there are ~ 1.6 times as many Compton-thick AGN as 
there are AGN with absorption column densities in the range 
10 23 < Nu < 10 24 crrT 2 . We apply the same correction here. 
Although very uncertain, the net effect is to increase the bolo- 
metric luminosity density by only a factor 1.48 (i.e. 0.17 dex). 
Varying the Compton-thick fraction from this nominal value 
would have an equivalent effect on the value of (e) deduced, 
but does not significantly affect the conclusions of this paper. 



4.2.3. The bolometric correction 

To convert 2 - lOkeV luminosity to bolometric lumin osity we 
use the bolometric correction of Ma rconi et alJ (2004), aver- 
aged over the distribution of luminosities of the llJeda et alJ 
(2003) model as a function of redshift. Because the luminosity 
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Fig. 5. The bol ometric lum i nosi ty density deduced from the 
best-fit model of llJeda et alJ d2003l) . integrating over the range 
10 40 < L x < lO^ergs -1 and applying the bolometric cor- 
rection of Marconi et al. (2004) and correction for Compton- 
thick AGN of lUeda et al.1 d2003l) (see text) (dashed curve). 
Also shown ar e uncertain t ies est imated from refitting to the 
binned data of lUeda et alJ (|2003) (points with error bars: see 
text; horizontal bars indicate the range of redshifts included 
in each point). For comparison, we show the e quivalent calcu- 
lation from the pure-luminosity-evolution fit of Richards et al. 
d2005l) to the optical 2SLAQ QSO survey in the redshift range 
0.4 < z < 2.1 (lower dashed curve). The luminosity density 
expected in PCE is shown for two cases: (i) no evolution in 
the comoving black hole mass density (dotted upper curve); 
(ii) evolution in the comoving black hole mass density that 
tracks the evolution of massive dark halos with Mh > 10 1I 5 M o 
(solid curve). Both curves assume average radiative efficiency 
(e) = 0.04 (see text). Note that the PCE calculation predicts a 
value for the luminosity density at z — 0, dependent only on 
(e) and cosmological parameters, with the higher redshift evo- 
lution further dependent only on choice of dark halo mass: the 
normalisation is not allowed to float arbitrarily. 



function is distributed to higher luminosities at higher redshifts, 
the average fraction of the bolometric luminosity that appears 
in the X-ray band varies from 0.056 at z = to 0.023 at z — 3, 
resulting in a variation in log I0 (pBOL) of 0.4 over this redshift 
range arising from the change in mean bolometric correction. 
Hence the bolometric luminosity density shows stronger evolu- 
tion, and has a maximum at higher redshift, than the 2 - lOkeV 
X-ray luminosity density alone. The bolometric correction is 
notoriously uncertain, and one key uncertainty is the non-linear 
dependence between X-ray and optical luminosity that leads 
to th e strong luminosity- dependence of the bolo metric correc- 
tion. [Marconi et all d2004l) assume the relation of lVignali et alJ 
( 200 3J), which has rece ntly been extended to lower luminosities 
bv lSteffen etafl d2006l) and which now seems fairly robust. To 
test this further we also compare with the bolometric luminos- 
ity density derived from optical QSO surveys. 
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4.2.4. The bolometric luminosity density from optical 
QSO surveys 

Optical QSO surveys do not probe to sufficiently low in- 
trinsic luminosities to pick up the significant population of 
Seyfert and type II AGN, a nd they are domina t ed by broad- 
line t ype I objects (e.g . iRichards etai1l2005t ISteffen et all 
2006 ; Bargeretal. 2005). However, it will be a useful test 
of the Marco ni et alJ IJ2004) bolometric correction to attempt 
to measure the bolometric luminosity density from an opti- 
cal QSO survey: the X-ray and optical bolometric corrections 
have a dependence on bolometric luminosity of opposite signs 
( Ma rconi et alll2004l) . arising from the L x - L opl correlation al- 
ready mentioned, so if the bolometric correction is wrong this 
will be manifest in this test. 



The optical QSO survey that to date has most success- 
fully probed_to_Jow intrinsic luminosity is the 2SLAQ sur- 
vey (Richard s et all 120051) . In Fig.|5]we plot the bolometric 
luminosity density obtained by integrating the pure luminos- 
ity evolution model derived from the 2SLAQ data alone, inte- 
grating over the equiva lent luminosity range as for the X-ray 
determination, with the Marconi e t alJ ij 2004) optical bolomet- 



ric correction and Uedaetal. (2003) Compton-thick fraction 
as above. Because of optical selection bias es the valid re dshift 
range is restricted to 0.4 < z < 2.1 ( Richard s et allE005t) . 

The function shown is more sensitive to the limits of in- 
tegration than the X-ray case above, because although the 
large survey areas allow optical QSO surveys to probe to 
higher AGN bolometric luminosities, the range of bolomet- 
ric luminosity covered is substantially smaller than the com - 
posite X-ray surveys discussed above (Richard s et all 12005). 
Increasing the lower luminosity limit to match a value Lx > 
10 42 ergs _1 decreases the z = 0.4 bolometric luminosity den- 
sity by 0.19dex and the maximum value by 0.14dex. Again, 
decreasing the upper luminosity limit to match a value Lx < 
10 46 ergs~' leads to no significant change in luminosity den- 
sity, to 0.01 dex. 

Considering the differing selection of AGN types in opti- 
cal and X-ray surveys, the X-ray- and optical-derived bolomet- 
ric luminosity densities show remarkably good agreement. This 
demonstrates that the corrections applied are unlikely to be too 
far wrong, at least over the optical-X-ray part of the spectrum 
(both bolometric corrections could be missing components in 
other parts of the s pectrum such as the far infrared - this is 
discussed further by M arconi et al.ll2004l) . The differences be- 
tween the two functions are most apparent at z ~ 0.7, although 
even here they agree within a factor 1.6. This difference arises 
from the additional low-luminosity AGN that are found by hard 
X-ra y surveys but not b y the large-area optical QSO surveys 
(see lRichards et al.l2005l for further discussion). 

This diagram also illustrates the limited effect of possi- 
ble "cosmic downsizing" on the luminosity density: the optical 
density is derived from a pure luminosity-evolution model, the 
X-ray from a luminosity-dependent density evolution model. 
We can see that the two models are measurably different, but 
that any "cosmic downsizing" does not dominate the redshift 
evolution of the luminosity density. 



4.2.5. The predicted bolometric luminosity density 

To calculate the luminosity density expected from PCE, we as- 
sume an average fraction (e) of the rest energy of accreted bary- 
onic material is radiated. Since the fractional accretion rate is 
almost independent of mass, we can write the expected bolo- 
metric luminosity density as 



<PL> 



c Pbh 



(1 



2 /e/(M H ) 
c Pbh 



e) M BH dt 

d5 c 



d-e) 



dt 



(19) 



where p B n is the cosmic black hole mas s dens ity. For the 
latter quantity we follow Mar coni et all J2004h and inte- 
grate t he black hole mass fu nction obtained by convolv- 
ing the INakamura et all J2003I) luminosity functio n with the 
luminosity-MBH relation of lMarconi & Hunl J2003I) . which in- 
tegrated for Mbh > 10 5 M o yields a z — mass density 
Pbh = 4.8 X 10 5 M o Mpc~ 3 with an uncertainty of about 30 
percent ( Marconi et all |2004 ) . 

The result is shown in Fig.[5]for two cases. The first shows 
what would happen if the black hole mass density did not 
change with cosmic epoch: the luminosity density would then 
be determined entirely by the invariant black hole mass den- 
sity and the epoch-dependent accretio n rate, as show n by the 
dotted curve. For consistency with llJeda et all J2003I) we as- 
sume a ACDM cosmology with parameters n« = 0.3, Hq = 
73km s _1 Mpc~ 3 (the observation al values are also ad justed to 
this value of Ho) an d °"8 = 0.74 ( Spe rgel et alll2006l) . The re- 
maining parameter is the average radiative efficiency (e) which 
we have set to a value 0.04 to approximately match data and 
prediction at z > 0.5. Such a value is consistent both with 
the range of values that might be expected for accretion onto 
black holes and with determinations that compare the X-ray 
background with the local black hole mass function, leading 
to a luminosity-function-dependent value for the radiative ef- 
ficiency specifically o f luminous AGN of (6agn) ~ 0.08^^ 
( Mar coni et all I2004T) . Note that this could differ from the av- 
erage efficiency for all black holes if a proportion were growing 
by radiatively-inefficient accretion: taken literally the agree- 
ment with the PCE model with (e) = 0.04 implies that such 
inefficient growth does not dominate the overall growth of 
black holes at z ^ 0.5, although some "ADAF" contribution 
is allowed. However, the deduced value of (e) is degenerate 
with the rather uncertain values of the bolometric correction 
and black hole mass density and with the value of <x 8 . The 
local black hole mass density is uncertain by ~ 30 percent 
and the value of <xs currently has an uncertainty ~ 10 percent 
(Sperg eTet all Eo()6). The uncertainty in the bolometric cor- 
rection is not well determined. Also, introduction of non-linear 
coupling between black hole and dark halo growth ("section RTTl 
would also modify the deduced value of (e) by a factor a (equa- 
tionll8l. Hence we should not attach too much importance to 
its value, provided that it is in the range expected for accretion 
onto supermassive black holes (e.g. ~ 0.06 for a Schwarzschild 
black hole). However, a significantly lower value might imply 
a significant contribution to black hole growth either from ra- 
diatively inefficient accretion or from obscured growth. 
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The second case shown in Fig.|3]considers also the expected 
evolution in the black hole mass function: if black holes do in- 
deed grow coevally with their host galaxies and halos we ex- 
pect their mass function to show similar evolution. The lower 
solid curve shows the result expected if the integrated black 
hole mass density evolves the same way as massive dark halos 
with Mh > 10 1l5 M<t). This e voluti on was calculated by inte- 
grating the ISheth & Tormenl J2002h mass function above this 
limit. 1 We find that imposing an upper limit has little effect, 
but that the amount of turn-down with increasing redshift in 
the integrated luminosity density does depend on the lower halo 
mass limit, with lower masses having less high-z turn-down, as 
expected from standard CDM models. 

4.2.6. Comparison of observed and predicted 
luminosity densities 

At z ^ 0.5 the simple hypothesis that black holes grow with 
their parent halos seems remarkably successful. With a fixed 
set of cosmological parameters, an independently determined 
estimate of the local mass density in black holes and a radiative 
efficiency in the expected range for AGN, this simple model 
of mass build-up predicts the local AGN luminosity density to 
within the observational uncertainties. No other tuning of the 
model or variation of parameters is required. 

Further, the integrated luminosity density is dominated by 
the contribution from AGN around the "break" in the luminos- 
ity function, with Lx ~ 10 44 ergs~': AGN of this l uminosity 
typic ally have black hole masses around 1O 8 M (IFine et all 
2006) and QSOs of this luminosity in the 2dF QSO Redshift 
Survey (2QZ) have dark halo masses independently estimated 
from their clustering bias of Mh ~ 10 125 /z _1 M (ICroom et all 
2005). Thus the higher-redshift evolution of the luminosity 
density is entirely in accord with the expectation that it is dom- 
inated by accretion onto black holes contained within massive 
galaxy halos. 

At z < 0.5 the luminosity density falls off faster than pre- 
dicted, although even at z = the predicted and observed val- 
ues agree to a factor two. The implication is that although halos 
and black holes grow coevally at higher redshifts, in the more 
nearby universe black hole growth may have decoupled from 
halo growth. There is growing evidence that accretion rate in 
the low-redshift universe depends on host galaxy mass, with 
lower mass galaxies having black holes that at low redshift are 
preferentiall y growing with respec t to their higher-mass coun- 
terparts ("e.g. lHeckman et a"fll2004l) - so-called "cosmic down- 
sizing". In detail, to explain cosmic downsizing we may need 
to recognise that in the late universe higher mass galaxies have 
slowed down their growth more than lower-mass galaxies, per- 
haps as a result of the environment they find themselves in at 
late epochs, or perhaps as a result of feedback. An alterna- 
tive explanation is that the mean radiative efficiency may fall 



if AGN accreting at lower Edding ton ratios switc h to an ADAF 
mode of accretion (e.g. lBeckert & Duschll l2002t) . 

What we can learn from the PCE model however is that 
these processes are secondary to the primary cause of AGN 
evolution: the cosmic build-up of structure and its rate of 
change. We reiterate that the success of the PCE hypothesis 
does not mean that complex physics is not operating: it sim- 
ply means that averaged across all massive galaxies, there is 
a mean net effect which corresponds to black holes growing at 
about the same rate as their host dark halos. Individually we ex- 
pect black holes to go through periods of inactivity punctuated 
by bursts of accretion probably associated with merger events 
between galaxies, remembering that the galaxy mass function 
itself is also built up hierarchically by mergers between less 
massive progenitors. But just as the average effect on the galaxy 
mass function is a steady build-up of mass, so there is a steady 
build-up of mass in the black hole population. 

A final consistency check is to ask whether the inferred 
local black hole mass density could indeed have been built 
up during the process of lu minous accretion. Th is question 
has largely be en answered by l|Marconi et alll20041) . who have 
shown that the lUeda et alJ J2003I) luminosity Sanction correctly 
predicts not only the hard extragalactic X-ray background spec- 
trum but also the local black hole mass density, as already dis- 
cussed above. The implication is that not only are the estimates 
of bolometric luminosity density consistent with the relic black 
hole population, but that there cannot have been any substan- 
tially larger amount of "hidden" (Compton-thick or radiatively- 
inefficient) black hole growth, unless the observed AGN have 
an uncomfortably high radiative efficiency (e » 0.1). 

4.3. The evolution in the Eddington ratio of 
supermassive black holes 

An alternative view of the evolution expected from the PCE hy- 
pothesis is to calculate the mean Eddington ratio {A), where A is 
the ratio of the actual mass accretion rate M acc to the Eddington 
mass accretion rate Mndd defined by LEdd = eMnddC 2 , so that 



(A) = 



/ eccr T M acc \ 
\4-nGM Bll mp I 
co-t I e/(M H ) 



4nGmp \ (1 — e) 



d5c 
dt 



(20) 



1 Note that it would be incorrect to attempt to calculate the evolution 
in Pbh(z) by integrating (a 'log M/dt): in hierarchical growth, mergers 
are important and evolution of the mass function is described by the 
Press-Schechter function. 



In deriving this expression, it is important to remember that 
we have calculated the PCE-predicted average Eddington rate 
of all black holes in galaxies, not simply those that happen to 
be visible as AGN at any particular epoch. As we don't expect 
accretion to be a smooth continuous process for any one galaxy, 
this is an important distinction. An individual galaxy with an 
observed high rate of mass accretion could have an individual 
Eddingto n ratio approac hing unity as seen in even local AGN 
fe.g. lOnken etalj l2004) but the whole population of galaxies 
would have a much lower mean Eddington ratio. 

So what significance can we attach to the mean Eddington 
ratio? If (A) has a high value, it implies that, on average, 
there is plenty of matter available to fuel luminous accreting 
black holes at or close to their Eddington rate. In fact, val- 
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Fig. 6. Ratio of mean accretion rate to the limiting Eddington 
rate, A, as a function of (1 + z) over the range < z < 20, 
plotted logarithmically to demonstrate the near power-law de- 
pendence. Results are shown for Mh = 10 12 5 h~ l M G and mean 
radiative efficiency (e) = 0.04 for: (a) ACDM cosmology, 
Q = 0.3, Q A = 0.7, H = 73kms-' Mpc~', cr 8 = 0.74 (solid 
curve); (b) Einstein-de Sitter cosmology with parameters as in 
(a) except £2q = 1, Qa = (dashed curve). 



ues significantly higher than unity would imply that there is 
actually too much accreting matter available: we should ex- 
pect that black hole growth by accretion would be limited to 
the Salpeter rate but that overall growth by mergers might be 
important. Conversely, a low value implies that, on average, 
there is insufficient new material accreting onto massive galax- 
ies to maintain Eddington-limited black hole growth. Fig. [6] 
shows (A) as a function of redshift for a flat, ACDM cosmol- 
ogy, assuming D.m - 0.3, Hq = 73kms~ 1 Mpc~', <x 8 = 0.74, 
M H = 10 I2 - 5 /i~' M , radiative efficiency e = 0.04 and without 
non-linear modification (section l4~Tt . 

The result plotted is remarkable in two respects. As already 
discussed, (A) depends only weakly on the mass of dark mat- 
ter halo and on the shape of the dark matter power spectrum, 
with no dependence on the fraction of a halo's total mass (dark 
matter plus baryons), that ends up in the black hole. Even the 
dependence on cosmology is weak, with a near power-law de- 
pendence on 1 + z with exponent n ^ 2.4 at z < 1 for the 
ACDM cosmology compared with n = 2.5 for Einstein-de 
Sitter. Second, in principle there was no reason why the growth 
of dark matter halos should have any interesting relationship to 
the Eddington limit, except that both are governed by gravity. 
Yet it turns out that at low redshifts the mean ratio (A) is con- 
sistent with the range of observationally determined values A L 
(determined as the ratio of bolometric luminosity to estimated 
Eddington luminosit y) for galaxy samples with a mean value 
(A L ) ~ 0.0035 aHcl |2004 as discussed below. At moderate 
redshifts the mean ratio increases to a value (A) ~ 0. 1 at z ~ 4, 
implying that at that epoch every galaxy has an actively accret- 
ing nuclear black hole (although typically with lower mass than 
today's black holes). 



5. Discussion 

5.1. Dispersion in Eddington ratio 

The above analysis has only enabled us to calculate the mean 
evolution in luminosity density and in (A). Since accretion 
events are discrete, we expect there to be a distribution of 
A values at any given cosmic epoch. In principle one could 

either attempt to calculate that distribution from EPS the- 

l 1 II I 

ory (e.g. Cole et al. 2000) or from numerical simulation (e.g. 
iHopkins et alJl2006llDi Matteo et alJE o05). or one coul d look 
at the observational constraints on the distribution. iHol (2004) 
has estimated the distribution of Al for Palomar galaxies, with 
bolometric luminosities estimate d from Ha and black hole 
mass estimated from applying the Tremai ne et alJ (|2002) rela- 
tion to measured stellar velocity dispersions. The distribution at 
essentially zero redshift is broad, with a mean {Al) ~ 0.0035, 
and although there are large uncertainties associated with the 
determination of Al the agreement to a factor about 2 with the 
value {A} ~ 0.0015 obtained at z _ for <e> = 0.04 (shown in 
Fig.|6) is encouraging. The lHol(l2004 values have a large dis- 
persion, which results in there being a smal l number of gal axies 
with A L ~ 0.1 that are visible as AGN (Onk en etall E004). and 
the vast majority of galaxies accreting at rates closer to (Al). 

We again emphasise that the mean value (A) calculated here 
is the mean for all galaxies and is not the mean value expected 
for samples of AGN. By definition AGN are selected because 
they are actively accreting and they must therefore be biased 
to higher values than the mean. Using the results from this pa- 
per to predict the distribution and evolution of AGN requires 
knowledge of or assumptions about the distribution of A and 
will be discussed further by Babic et al. (in preparation). 

5.2. AGN evolution and accretion at high redshift 

As seen in section l4~2l the absolute value and the evolution 
of either the observed space or luminosity density appears to 
be a combination of two competing but inter-related effects: 
the build-up of massive halos and hence black holes with de- 
creasing redshift, and the decreasing accretio n rate with de- 
creasing redshift calculated in this paper. Icroom et all |2005) 
have shown that luminous QSOs inhabit dark matter halos of 
mass ~ 10 12 5 /i _1 M : at low redshifts (z < 1) the halo mass 
function, and hence the black hole mass function, is chang- 
ing little at this mass, and AGN evolution is dominated by 
the evolution in A. In this case we expect the evolution of 
the population to be dominated by apparent luminosity evo- 
lution: on average black holes at lower redshift accrete at a 
lower rate than black holes a t higher redshift. B roadly speak- 
ing, this is what is observed ( Bovle et al. , 2000), although it's 
not the whole sto ry (Steffe n et all 120031 ICowie et all 120031 
lUeda et all Eoolt IZheng et all l2004 iBarger et all 120051) . At 
higher redshifts (z > 2) the growth of the black hole mass func- 
tion is significant and probably dominates the evolution, caus- 
ing a decline in space density with redshift (Fig.|5Jl. 

At z ~ 4 (A) approaches the values found for the most lu- 
minous QSOs at lower z: at this epoch in cosmic history the 
average dark matter halo is accreting at a rate sufficient to sup- 
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ply close to the Eddington luminosity of its nuclear black hole: 
the majority of galaxies with such a black hole would be ac- 
tive. At higher redshifts still the rate of accretion of matter onto 
dark halos exceeds that required to supply the Eddington rate. 
In this case it is unlikely that super-Eddington accretion onto 
black holes would occur: it is more likely that nuclear outflows 
would limit the accretion process of individual black holes to 
about the Eddington rate and it is also possible that the black- 
hole/bulge relation may be built up during this phase of ga laxy 
and black hole growth (King, 2003:lKin g & Pounds! Eotm 



5.3. Comparison with previous work 
In previous numerical (e.g. iKauffmann & Haehneltl 1200 



Volonteri et all 120031 ICroton et all l200fj) and analytic (e.g. 
Percival & Miller! 1 1999) attempts to understand AGN evolu- 



tion and supermassive black hole growth, AGN activity was 
hypothesised to be directly associated with merger events, and 
prescriptions were adopted that describe the relationship be- 
tween the halo mass of a galaxy and the resulting AGN's lumi- 
nosity, and between the rate of mergers and the evolving space 
density of AGN. Implicitly, a relationship between black hole 
mass and halo mass is assumed. To obtain the correct space 
density of AGN it was also necessary to assume a timescale for 
the luminous phase of a black hole following a merger. In these 
models, the evolution in EPS merger rate effectively produces 
number density evolution, whereas the dominant component 
of observed QSO evolution takes the form of luminosity evo- 
lution. As argued in the introduction, there must in fact have 
been evolution in mean accretion rate if local massive galaxies 
contain dormant supermassive black holes. This work shows 
how that evolution arises naturally in ACDM models. There 
has been no requirement to make assumptions about the lumi- 
nosity/halo relation or the link with mergers and no assump- 
tion required about the timescale of AGN activity. The calcula- 
tion presented here has assumed coeval growth of bl ack holes: 
similar assumptions are made in, for example, Volonteri et alJ 
J2003I) . but this assumption certainly needs further testing and 
evaluation. Overall, the extremely good agreement between the 
predicted and observed AGN luminosity density and its evo- 
lution, without any "model-tuning", points to a close link be- 
tween AGN cosmological evolution and the evolution in accre- 
tion onto galaxy halos. Conversely, predictions of recent semi- 
analytic models that attempt to provid e recipes for the com plex 
physics governing black hole growth (Crot on et al.!l2006l) pre- 
dict a peak in black hole formation rate (and hence presumably 
emitted luminosity density arising from accretion) at z ~ 3, 
with a decline to z — of a factor ~ 6, in disagreement with the 
observed luminosity density presented in section l4~2l 

Numerical simulations have also been used to try to un- 
derstand the link between black hole growth and galaxy merg- 
ers, and in particular to understand the role of feedback in 
forming the well-defined relation be tween galaxy velocity 
dispersion and black hole mass (e.g. | R obertso n et alJ |2006; 
IPi Matteo et all 12005! iHopkins et all 12005! l2006!) ~These ap- 
proaches are directly complementary to the question addressed 
in this paper, namely to understand analytically the effect that 



the cosmic slowdown in halo growth has on the typical black 
hole accretion rate. It may be that combining the two ap- 
proaches will finally lead us to a detailed understanding of 
AGN evolution. 



6. Conclusions 

The results from this work may be summarised as follows. 

1. Press-Schechter theory may be used to calculate analyti- 
cally the mean mass accretion rate of dark matter halos. 
The result of the calculation is in good agreement with 
previous numerical work but offers a significant improve- 
ment over previously-published fits to numerical simula- 
tions. The mass accretion rate depends almost linearly on 
mass but is very insensitive to choice of cosmological pa- 
rameters within the range normally considered. We present 
fitting formulae to allow easy calculation of the accretion 
rate. 

2. It seems likely that the deduced strong decline in halo mass 
accretion will have an effect not only on the build-up of 
dark halos but also on the baryonic structures they contain. 
There may therefore be a close link between the cosmolog- 
ical evolution in star formation rate, in AGN accretion and 
dark halo accretion. 

3. We have investigated arguably the simplest hypothesis one 
could make for black hole growth, "Pure Coeval Evolution" 
(PCE), which postulates that, on average, supermassive 
black holes have growth that tracks the growth of their par- 
ent dark halos. In this case it is straightforward to calcu- 
late the expected mass accretion rate onto black holes and, 
estimating the local black hole mass density, we predict a 
value for the integrated AGN luminosity density in remark- 
able agreement at z ^ 0.5 with the value deduced from hard 
X-ray surveys. 

4. The evolution in integrated AGN luminosity density at 
z > 0.5 is well-matched by this model if the integrated 
black hole mass density tracks the mass density in massive 
Mh > 10 115 M G halos, consistent with the value obtained 
by inferring th e clustering b i as of QSOs in the 2dF QSO 
redshift survey JCroometall E005) 

5. At z < 0.5 the observed luminosity density falls off faster 
than predicted, although even by z — the overprediction 
is only a factor two. It seems that in the recent universe the 
black hole growth has started to decouple from dark halo 
growth, although a decrease in mean radiative efficiency as- 
sociated with an increasing prevalence of ADAFs may also 
occur. 

6. Expressed in terms of the Eddington ratio, at sufficiently 
high redshifts (z S 4), halos would have black holes 
accreting at a significant fraction of the Eddington rate 
- we expect all massive galaxies to contain an active 
AGN. At higher redshifts still, the nominal Eddington ratio 
could greatly exceed unity. The response of the black hole 
would be to limit its growth to the Salpeter r ate, but such 
systems might produc e significant outflows llKing! 12003! 
King & Pounds, 2003) which in turn would produce a sig- 
nificant feedback effect on their hosts. It is also likely that 
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in this phase of black hole growth, growth by mergers be- 
tween halos is faster than growth by accretion. Merger- 
dominated growth is one way of avoiding the problem of 
high-mass black holes at high redshift that have had in- 
sufficient time to grow at the Salpeter rat e from solar- 
or int ermediate-mass black hole progenitors JWillott et all 
120031) . 

Further work will model the AGN luminosity function and the 
contribution to the X-ray background (Babic et al. in prepara- 
tion). 
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Appendix A: Calculating dd c /dt 

In this section we discuss the calculation of dS c /dt and approxima- 
tions to it for cosmologies with a cosmological constant (the analysis 
could be extended to include more general dark energy cosmologies 
following the approach of Percival (2005) but this is beyond the scope 
of this appendix). We consider quantities that evolve as a function of 
scale factor a = (1 + z)" 1 . If no dependence is quoted, the quantities 
should be assumed to be calculated at a = 1 (z = 0). 

The critical overdensity for collapse S c (a) in the spherical top-hat 
collapse model is the linear overdensity, extrapolated to the present 
day, that leads to the collapse of a homogeneous spherical region to a 
singularity at scale factor a. This is used in EPS theory to link over- 
densities with their predicted collapse times. However, two alternative 
formalisms are often considered for the collapse of perturbations in 
EPS theory. 

1. The overdensity field is assumed to grow with the linear growth 
factor, and when perturbations reach a particular critical over- 
density they are said to have collapsed. This suggests that, given 
5 C , we can calculate S c (a) = 6 C D/D(a), where D(a) is the linear 
growth factor. 

2. An overdense region is considered to be spherical and to evolve 
according to the Friedmann equations, and hence its collapse time 



14 



L. Miller et al.: Cosmological accretion history 




Fig. A.l. Plot showing the mismatch between contours of con- 
stant 8 C (black lines) and the lines along which Q.m and Qa vary 
with the evolution of the universe (grey lines). The contours of 
constant 5 C are plotted as a ratio of c>Eds, the critical overdensity 
for collapse at redshift zero in an Einstein-de Sitter cosmology. 
Because these lines cross, the evolution of 6 c (a) is not solely 
dependent on the linear growth factor. 
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Fig.A.2. The percentage error in approximating the time 
derivative of 6 c (a) as 6 c (a) = 6EdsD/D{a), where D(a) is the 
linear growth factor. Open symbols incorporate the approxima- 
tion for D(a) given by Eq. IA.6I while solid symbols calculate 
D(a) by numerically integrating Eq. IA.4I Results are shown for 
three cosmological models: for an Einstein-de Sitter cosmol- 
ogy, both approximations reduce to the exact behaviour. 



may be calculated. The overdensity value is extrapolated to the 
epoch at which the density field is normalised (conventionally, the 
present day) using the linear growth factor. 

These two alternatives give different formulae for 6 c (a). In 
Fig. IA.ll we plot contours of constant S c (black lines) w ithin the p lane 
of (A m,Qa) va l ues > calculated using the method of IPercival et alJ 
(2000). Contours are plotted as a function of the value of S c in an 
Einstein-de Sitter cosmology: 

<5 Ed s = ^(12^) 2/3 - 1.686. (A.l) 
We also plot lines along which universes evolve (grey lines) given by 



Q M (a) 



£2 A (a) 



a + (l 



a)Q.M + (a 3 
a 3 Q A 



(A.2) 



(A.3) 



a + (1 — a)Q.M + {a — a)Q. A 

Because 6 C changes along the lines of evolution, we see that the spher- 
ical top-hat collapse model does not predict evolution of the critical 
overdensity given by S c (a) = 6 C D/D(a). 

The first picture of EPS theory as corresponding to a grow- 
ing field does not follow from the spherical top-hat collapse model. 
However, Fig. lA. ll shows that the evolution of S c along lines of evolv- 
ing (fl M (a), Sl A (a)) is small. In the remainder of this appendix we con- 
sider the error introduced in \dS c (a)/dt\ by assuming that the evolu- 
tion in 6 c (a) is governed solely by the linear growth factor. The linear 
growth factor D(a) is given by 



D(a) = -a M HlH(, 



a) f 

Jo 



da 



(aH(a)f 



(A.4) 



where 

<l_da\ 
\a dt 



= H(af = Hl(a K + n M cc 3 - (Q,„ + n A - ik 2 ). (a.5) 



D(a) may be approximated using the fitting formula of Carroll et al. 
(1992): 
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Eq. lA.4l can be differentiated to obtain 
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Finally, the approximation to dS c (a)/dt is given by 



dSAfl) 
dt 



D dD(a) da 
; £»(a) 2 da It' 



(A.7) 



(A.8) 



Eauations lA.5IIA.6llA,7l & lA,8l combine to provide an analytic ap- 
proximation to dS c ldt. In Fig. lA. 21 we plot the percentage error intro- 
duced by this approximation compared with numerically integrating 
the behaviour of S c (a), calculated using the method of Pe rcival et alJ 
12000) (open symbols). We also consider calculating D(a) exactly 
and use this instead of the Carroll et al. (1992) approximation in 
Ea. lA.7l (solid symbols). The solid symbols show the error in dd c /dt 
from ignoring the evolution in S c because of changing C1 M an d ^a- 
Comparison of solid and open symbols shows the error introduced by 
Eq. IA.6l Errors are plotted for three cosmologies: for the Einstein-de 
Sitter cosmology, the approximations reduce to the exact result. For 
the A and open cosmologies, the maximum error is less than 2%. 



